An Extended Scheme for Fitting X-ray Data with Accretion Disk 
Spectra in the Strong Gravity Regime^ 
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ABSTRACT 

Accreting black holes arc believed to emit X-rays which then mediate information about strong 
gravity in the vicinity of the emission region. We report on a set of new routines for the xspec 
package for analysing X-ray spectra of black-hole accretion disks. The new computational tool 
significantly extends the capabilities of the currently available fitting procedures that include the 
effects of strong gravity, and allows one to systematically explore the constraints on more model 
parameters than previously possible (for example black-hole angular momentum). Moreover, 
axial symmetry of the disk intrinsic emissivity is not assumed, although it can be imposed to 
speed up the computations. The new routines can be used also as a stand-alone and flexible 
code with the capability of handling time-resolved spectra in the regime of strong gravity. We 
have used the new code to analyse the mean X-ray spectrum from the long XMM-Newton 2001 
campaign of the Seyfert 1 galaxy MCG-6-30-15. Consistent with previous findings, we obtained a 
good fit to the broad Fe K line profile for a radial line intrinsic emissivity law in the disk which is 
not a simple power law, and for near maximal value of black hole angular momentum. However, 
equally good fits can be obtained also for small values of the black hole angular momentum. 
The code has been developed with the aim of allowing precise modelling of relativistic effects. 
Although we find that current data cannot constrain the parameters of black-hole/accretion disk 
system well, the code allows, for a given source or situation, detailed investigations of what 
features of the data future studies should be focused on in order to achieve the goal of uniquely 
isolating the parameters of such systems. 



Subject headings: black hole physics — line: profiles - 
individual (MCG-6-30-15) — methods: data analysis 

1. Introduction 

There is now strong evidence that the Fe Ka 
line emission in some active galactic nuclei (AGN) 
and some Galactic X-ray binary black-hole can- 
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X-rays: galaxies — galaxies: active — galaxies: 

didates (BHC) originates, at least in part, from 
an accretion disk in a strong gravitational field. 
A lively debate is aimed at addressing the ques- 
tion of what the spectral line profiles and the as- 
sociated continuum can tell us about the central 
black-hole, and whether they can be used to con- 
strain parameters of the accretion disk in a nearby 
zone, about ten gravitational radii or less from 
the center. For a recent review for AGNs, see 
Fabian et al. (2000), Reynolds & Nowak (2003), 
and references cited therein. For BHCs, see Miller 
et al. (2002a), McClintock & Remillard (2003) 
and references therein. In several sources there 
is indication of Fe Ka line emission from within 
the last stable orbit of a Schwarzschild black hole 
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(e.g. in the Seyfert galaxy MCG-6-30-15; see Iwa- 
sawa et al. 1996; Fabian et al. 2000; Wilms et al. 
2001; Martocchia, Matt & Karas 2002a) while in 
other cases the emission appears to arise farther 
from the black hole (e.g. in the microquasar GRS 
1915+105, sec Martocchia et al. 2002b; for AGNs, 
see Yaqoob & Padmanabhan 2004 and references 
therein). Often, the results from X-ray line spec- 
troscopy are inconclusive, especially in the case of 
low spectral resolution data. For example, a spin- 
ning black hole is allowed but not required by the 
line model of the microquasar V4641 Sgr (Miller 
et al. 2002b). The debate still remains open, but 
there are good prospects for future X-ray astron- 
omy missions to be able to use the Fe Ka line to 
probe the space-time in the vicinity of a black hole, 
and in particular to measure the angular momen- 
tum, or the spin, associated with the metric. 

One may also be able to study the 'plunge 
region' (about which very little is known), be- 
tween the event horizon and the last stable or- 
bit, and to determine if any appreciable contri- 
bution to the Fe Ka line emission originates from 
there (Reynolds & Begelman 1997; Krolik & Haw- 
ley 2002). In addition to the Fe K lines, there is 
some evidence for relativistic soft X-ray emission 
lines due to the Lya transitions of oxygen, nitro- 
gen, and carbon (e.g. Mason et al. 2003), although 
the observational support for this interpretation is 
still controversial (e.g. Lee et al. 2001). 

With the greatly enhanced spectral resolution 
and throughput of future X-ray astronomy mis- 
sions, the need arises for realistic theoretical mod- 
els of the disk emission and computational tools 
that are powerful enough to deal with complex 
models and to allow actual fitting of theoretical 
models to observational data. It is worth noting 
that some of the current data have been used to 
address the issue of distinguishing between differ- 
ent space-time metrics around a black hole, how- 
ever, the current models available for fitting X-ray 
data are subject to various restrictions which we 
shall elaborate on in the present paper. 

In this paper we describe a generalised scheme 
and a code which can be used with the stan- 
dard X-ray spectral fitting package, xspec (Ar- 
naud 1996). We have in mind general relativity 
models for black- hole accretion disks. Apart from 
a better numerical resolution, the principal inno- 
vations compared to currently available schemes 



are that the new model allows one to (i) fit for 
the black-hole spin, (ii) study the emission from 
the plunge region, and (iii) specify a more general 
form of emissivity function of polar coordi- 
nates in the disk plane (both for the line and for 
the continuum). Furthermore, it is also possible to 
(iv) study time variability of the observed signal 
and (v) compute Stokes parameters of a polarized 
signal. Items (i)-(iii) are immediately applicable 
to current data and modelling, while the last two 
mentioned features are still mainly of theoretical 
interest at present. Time-resolved analysis and po- 
larimctry of accretion disks are directed towards 
future applications when the necessary resolution 
and the ability to do polarimctry are available in 
X-rays. Thus our code has the advantage that it 
can be used with time-resolved data for reverber- 
ation studies of relativistic accretion disks (Stella 
1990; Reynolds et al. 1999; Ruszkowski 2000; Goy- 
der & Lasenby 2004). Also polarimetric analysis 
can be performed, and this will be extremely use- 
ful because it can add very specific information on 
strong-gravitational field effects (Connors, Stark 
& Piran 1980; Matt, Fabian & Ross 1993; Bao 
et al. 1997). Theoretical spectra with temporal 
and polarimetric information can be analysed with 
the current version of our code and such analysis 
should provide tighter constrains on future models 
than is currently possible. 

As mentioned above, Fe K lines have been re- 
ported in the X-ray spectra of numerous AGNs 
and BHCs. These sources frequently exhibit re- 
markable variability patterns which are still diffi- 
cult to understand. Among these puzzling objects, 
the Seyfert 1 galaxy MCG-6-30-15 is radically dis- 
tinctive with its broad and skewed Fe Ka feature 
that persists in observations in spite of substantial 
variablity of the continuum. Interpretation of the 
Fe K line in terms of reflection from a relativis- 
tic black-hole disk has been found to be rather 
robust, but it has not been possible to fully ex- 
plore the model parameter space. In particular, 
constraints on a have only been obtained by fit- 
ting 'inverted photon data' which is made from the 
real data in a model-dependent way that already 
assumes a certain value of a. The puzzling stabil- 
ity of this line has been attributed, amongst other 
things, to general relativity effects very close to the 
black hole (Miniutti et al. 2003; Fabian & Vaughan 
2003). It is thus interesting to know if and how 
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a generalisation of the standard disk-line scheme 
together with accurate computation of general rel- 
ativistic effects can add new pieces of information 
or whether a more substantial modification of the 
whole picture is needed (for further discussion and 
references, see e.g. Weaver & Yaqoob 1998; Kro- 
lik 1999; Hartnoll & Blackman 2001; Dumont ct 
al. 2002). Different kinds of such generalisations 
have been proposed. Previous papers indicated 
that details of the model do matter. However, 
a strong gravitational field has almost invariably 
been required. In the present paper, we adopt the 
assumption that the line and continuum emission 
are produced by an irradiated disk near a rotating 
black hole, we relax some of the previous assump- 
tions, and we search for best-fitting parameters of 
the model. 

In the § 2 we describe the main features of the 
new computational tool. We list several variants 
of the code which address different problems of 
fitting X-ray spectra of black hole plus accretion 
disk systems. In § 3 we use the code to anal- 
yse time-averaged data for MCG-6-30-15 XMM- 
Newton CCD data as detailed in Fabian et al. 
2002). We point out that more complex models 
may not require a large value of the black-hole an- 
gular momentum, although the current data are 
consistent with a maximally rotating black hole. 
Finally in § 4, we summarise our results and state 
our conclusions. 

2. Computational Technique 

The new model, ky, is suited for use with the 
xspec package (Arnaud 1996). Several mutations 
of ky were developed, with an eye to specifically 
provide different applications, and linked with a 
common ray-tracing subroutine, which therefore 
does not have to be touched when the intrinsic 
(local) emissivity function in the model is changed. 

When the relativistic line distortions arc com- 
puted, the new model is more accurate than the 
laor model (Laor 1991) and faster than kerrspec 
(Martocchia, Karas & Matt 2000). These are the 
other two xspec models with a similar usage (see 
also Pariev & Bromley 1998; Gierlinski, Maciolek- 
Niedzwiecki & Ebisawa 2001; Schnittman & 
Bertschinger 2003). It is also important to com- 
pare the results from fully independent relativistic 
codes since the calculations are sufficiently corn- 



Table 2: Basic versions of the model. 



Name Typet Usage 

KYGIine additive Relativistic spectral line 

from a black hole disk. 

KYHrefl additive Compton reflection with an 

incident power-law (or a 
broken power-law) contin- 
uum. 

KYLcr additive Lamp-post Compton reflec- 

tion model. 

KYSpot To be used as Time-dependent spectrum 
an independent of a pattern co-orbiting 
code outside with the disk. 

XSPEC 

KYConv convolution Convolution of the rela- 

tivistic kernel with intrinsic 
emissivity across the disk. 

t Different model types correspond to XSPEC syntax and 
are defined by the way they act in the overall model and 
form the final spectrum. According to the usual conven- 
tion in XSPEC, additive models represent individual emis- 
sion spectral components which may originate e.g. in dif- 
ferent regions of the source. Additive models are simply 
superposed in the total signal. Multiplicative components 
(e.g. hrefl, discussed in the Appendix B) multiply the cur- 
rent model by an energy-dependent factor. Convolution 
models modify the model in a more non-trivial manner. 
See Arnaud (1996) for details. 



plex that significant differences can arise. Our 
code is more general than the currently available 
alternatives in XSPEC, so it can be used to assess 
the shortcomings of the available codes. We com- 
pared calculated profiles of the relativistic line in 
our model with several models available in litera- 
ture. In particular, our time-averaged profiles are 
in perfect agreement with those shown in Bcckwith 
& Done (preprint, 2004). Small (but potentially 
significant) distinctions from other models can be 
explained by omission of light-focusing effect in 
diskline and by insufficient resolution in the laor 
model. The discrepancies become more prominent 
in extreme situations of large inclination and/or 
very narrow ring. Furthermore, as discussed in 
the text below, various uncertainties in the form 
of intrinsic emissivity and especially in the dark- 
ening law introduce changes in spectrum which 
could not be examined by older models. For more 
comparisons see Appendix A. 

As far as the continuum component is con- 
cerned, an approximation was adopted in order 
to account for Compton reflection. To this end, 
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Table 1: Basic features of the new model in comparison with other black hole disk-line models. 



Model 


Effects that are taken into account 




Reference 


Energy shift of pho- 
tons/Lcnsing effect 


Black hole 

angular 

momentum 


Axisymmetry 
is assumed 


Steady source 
is assumed 


diskline 


yes/no 





yes 


yes 


Fabian et al. (1989) 


laor 


yes/yes 


0.998 


yes 


yes 


Laor (1991) 


kerrspec 


yes / yes 


<o,i) f 


no 


yes 


Martocchia et al. (2000) 


ky 


yes/yes 


<o,i> 


no* 


no 


This paper 



t The value of dimension-less a parameter is kept frozen. 

t A one-dimensional version is available for the case of an axisymmetric disk. In this axisymmetric mode, ky still 
allows a and other relevant parameters to be fitted (in which case the computational speed of ky is then comparable 
to laor). The results can be more accurate than those obtained with other routines because of the ability to tune 
the grid resolution. 



we employ a more basic approach than the pexrav 
model (Magdziarz & Zdziarski 1995). Simplic- 
ity of the adopted scheme affords a much higher 
computation speed than other Compton reflection 
models in xspec. Further details of the method 
and the calculations for the ray-tracing are given 
in Appendix A, and details of the Compton reflec- 
tion model we use are given in Appendix B (see 
also Dovciak et al. 2004). 

Among its useful features, the ky model allows 
one to fit various parameters such as black-hole 
angular momentum (a), observer inclination an- 
gle relative to the disk axis (9 ), and the size and 
shape of the emission area on the disk, which can 
be non-axisymmetric (see Table 1 for a summary 
of the basic features). A straightforward modifi- 
cation of a single subroutine suffices to alter the 
prescription for the disk emissivity, which is speci- 
fied cither by an analytical formula or in a tabular 
form. Our code allows one to change the mesh 
spacing and resolution for the (two-dimensional) 
polar grid that covers the disk plane, as well as 
the energy vector (the output resolution is even- 
tually determined by the detector in use when the 
model is folded through the instrument response) . 
Hence, there is sufficient control of the (improved) 
accuracy and computational speed. 

Furthermore, ky can be run as a stand-alone 
program (detached from xspec). In this mode 
there is an option for time-variable sources such 
as orbiting spots, spiral waves or evolving flares 
(e.g. Czerny et al. 2004, who applied a similar ap- 



proach to compute the predicted rms variability in 
a specific flare/spot model). The improved accu- 
racy of the new model has been achieved in several 
ways: (i) photon rays are integrated in Kerr in- 
going coordinates which follow principal photons, 
(ii) simultaneous integration of the geodesic devi- 
ation equations ensures accurate evaluation of the 
lensing effect, and (iii) non-uniform and rather fine 
grids have been carefully selected. 

A new and efficient code has been desirable in 
order to be able to perform spectral fits with suffi- 
cient accuracy and to deal with non-axisymmetric 
geometry of accretion flows. Also, obscuration ef- 
fects along the line-of-sight need to be taken into 
account, as well as the effects of relativity which 
act on photons, as they propagate through curved 
spacetime towards a distant observer. A reason- 
able trade-off must be achieved between mathe- 
matically elegant approaches (as initiated by de 
Felice, Nobili & Calvani 1974, and Cunningham 
1975) and straightforward numerical ray-tracing 
of a sufficiently huge number of photons (Karas, 
Vokrouhlicky & Polnarev 1992; Bromley, Chen & 
Miller 1997). The former method is less flexible as 
far as the source geometry is concerned. On the 
other hand, the latter method is computationally 
demanding (this problem can be partly solved by 
analytical integration of the rays in terms of el- 
liptical integrals, which is amenable in the Kerr 
metric, but only at cost of substantial complexity 
of the resulting code; cf. Martocchia et al. 2000; 
Cadez et al. 2003). Various strategies and ap- 
proximations have been developed by different au- 
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thors (Asaoka 1989; Laor 1991; Kojima 1991; Bao, 
Hadrava & 0stgaard 1994; Dabrowski & Lascnby 
2001; Fanton et al. 1997; Matt et al. 1993a, b; 
Semerak, Karas & de Felice 1999; Viergutz 1991; 
Zakharov 1994). Taking into account the practical 
experience with various approaches, we optimised 
the new routines with respect to rather antago- 
nistic requirements on speed, accuracy, generality 
and flexibility. 

Several versions of the routine have been pre- 
fabricated for different types of sources (Table 2): 

(i) an intrinsically narrow line produced by a disk, 

(ii) a relativistically blurred Compton-reflection 
continuum including a primary power-law compo- 
nent, (iii) the lamp-post model emissivity (Mar- 
tocchia et al. 2000), and (iv) time-dependent spec- 
trum of an orbiting or a free-falling spot. Default 
parameter values for the line model correspond to 
those in the laor model, but numerous options have 
been added. For example, in the new model one 
is able to set the emission inner radius below the 
marginally stable orbit, r- m < r ms (a). One can 
also allow a to vary independently, in which case 
the horizon radius, rh(a) = 1 + \/l — a 2 , has to 
be, and indeed is, updated at each step of the fit 
procedure. We thus define emission radii in terms 
of their offset from the horizon. Several arguments 
have been advocated in favour of having ri n ^r ms 
for the disk emission, but this possibility has never 
been tested rigorously against observational data. 
The set of ky-routines introduced above provide 
the tools to explore black-hole disk models and to 
actually fit for their key parameters, namely, a, 
9 , and r; n . Let us remind the reader that, apart 
from the model parameters which are subject to 
fitting, one may also specify the functional form 
of intrinsic emissivity of the disk as a function of 
polar coordinates and photon energy. 

As an example, an time-independent intrin- 
sic emissivity Iioc(E\ oc ,r) cx J? bc r ° r" a = was con- 
volved with our ray-tracing routine to produce a 
power-law type continuum component with rela- 
tivistic effects. Compton reflection was taken into 
account for intrinsic emissivity in the same ap- 
proximation as adopted in the (non-relativistic) 
model href I. The spectrum of the resulting rela- 
tivistic model, kyhrefl, therefore resembles the out- 
come of a multiplicative model, hrefl*powerlaw in 
xspec, which has been blurred by the relativis- 
tic disk kernel. It is worth noting that the effect 



of general relativity is practically negligible in the 
case of a locally power-law continuum because of 
the overall smearing of the spectrum across the 
disk surface. A non-trivial part in kyhrefl there- 
fore concerns the Compton reflection component 
and its mutual normalization with respect to the 
primary powerlaw originating at different radii in 
the disk. 

We have also produced a convolution-type 
model, kyconv, which can be applied to any ex- 
isting XSPEC model of intrinsic X-ray emission 
(naturally, a meaningful combination of the mod- 
els is the responsibility of the user). We remind 
the reader that kyconv is substantially more pow- 
erful than the usual convolution models in XSPEC, 
which are defined in terms of one-dimensional 
integration over energy bins. Despite the fact 
that kyconv still uses the standard xspec syn- 
tax in evaluating the observed spectrum (e.g. ky- 
conv(powerlaw)), our code performs a more com- 
plex operation. It still performs ray-tracing across 
the disk surface, so the intrinsic model contribu- 
tions are integrated from different radii. The price 
that one has to pay for the enhanced functionality 
is a higher demand for computational power. 

As mentioned above, time-dependent spectra 
are expected to eliminate a great deal of uncer- 
tainty in studying accretion disks (e.g. Done & 
Gierlihski 2003). As an example of using ky out- 
side of XSPEC, Fig. 1 represents an evolving spec- 
tral line in a dynamical diagram: a localized spot, 
free-falling from the marginally stable orbit down 
towards the black-hole horizon. Intrinsic emissiv- 
ity is assumed to fade gradually with the distance 
from the spot centre, with characteristic diameter 
of the spot - OAGM/c 2 {[3 = 100 in kyspot), but 
it is kept constant in time. Therefore the figure 
reflects the effect of changing energy shift along 
the spot trajectory. The in-spiral motion starts 
just below the marginally stable orbit and it pro- 
ceeds down to the horizon, maintaining the specific 
energy and angular momentum of its initial circu- 
lar orbit. One can easily recognise that variations 
in rcdshift could hardly be recovered from data 
if only time-averaged spectrum were available. In 
particular, the motion above r ms (a) near a rapidly 
spinning black hole is difficult to distinguish from 
a descent below r ms (0) in the non-rotating case. 

Before time-resolved spectra are available with 
sufficient quality, preliminary considerations can 
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Fig. 1. — Dynamical profile of a line produced by an in-spiralling spot in free fall. Left: Trajectory (in the 
equatorial plane) of the spot, defined by a constant angular momentum during the infall. Right: Spectrum 
of the spot. Energy is on the ordinate, time on the abscissa. The horizontal range spans the interval of 
2.4 orbital periods at the corresponding initial radius, r = 0.93 r ms (a). Here, the dimension-less angular 
momentum parameter of the black hole is a — 0.9, observer inclination is 8 = 45°. The observed photon 
flux is colour-coded (logarithmic scale with arbitrary units). Gradual decay of the signal and an increasing 
centroid redshift can be observed as the spot completes over two full revolutions and eventually plunges into 
the black hole. 



be based on the extremal values, and the range, of 
energy shift, g + and g_ , as well as the relative time 
delay which photons experience when arriving at 
the observer's location from different parts of the 
disk. This is particularly relevant for some narrow 
lines whose redshift can be determined more accu- 
rately than if the line is broad (e.g. Turner et al. 
2002; Guainazzi 2003; Yaqoob et al. 2003). Care- 
ful discussion of g± can, in principle, circumvent 
the uncertainties which are introduced by uncer- 
tain form of the intrinsic emissivity of the disk and 
yet still constrain some of parameters. Advantages 
of this technique were pointed out already by Cun- 
ningham (1975) and it was further developed by 
Pariev, Bromley & Miller (2001). 

Fig. 2 shows the extremal values of the red- 
shift factor for 9 = 30° (which often seems to 
be preferred by the observational data). These 
were computed along r = constant circles in the 
disk plane, together with the corresponding con- 
tour lines of g± = constant in the plane a versus 
r. In other words, radiation is supposed to origi- 
nate from radius r in the disk, but it experiences 
a different redshift depending on the polar angle. 
Contours of the redshift factor provide a very use- 



ful and straightforward technique to determine the 
position of a flare or a spot, provided that a narrow 
spectral line is produced and measured with suf- 
ficient accuracy. The most direct use of extremal 
<?-values would be, if one were able to measure 
variations of the line profile from its lowest-energy 
excursion (for gj) and highest-energy excursion 
(for g+), i.e. over a complete cycle. Even partial 
information can help to constrain models (for ex- 
ample, the count rate is expected to dominate the 
observed line at the time that g = g+, and so 
the high-energy peak is easier to detect). While 
it is very difficult to achieve sufficient precision 
on the highly shifted and damped red wing of a 
broad line, prospects for using narrow lines is in- 
deed interesting, provided that they originate close 
enough to the black hole and that sufficient res- 
olution is achieved both in the energy and time 
domains. 

Fig. 2 also makes it clear that it is possible in 
principle (but intricate in practice) to deduce the 
a-parameter value from spectra. It can be seen 
from the redshift factor that the dependence on a 
is rather small and it quickly becomes negligible if 
the light from the source is dominated by contri- 
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Fig. 2. — Energy shift of photons originating from different radii in the disk. The dependence on the black- 
hole angular momentum, a, is shown. Left: extremal values of the redshift factor g±(a\ r), (i.e. the maximum 
and minimum values, corresponding to two bundles of neighbouring curves). Each curve corresponds to a 
fixed value of a in the range (0, 1), encoded by line colours. Right: Corresponding contour lines of extremal 
energy shift factors in the plane of the rotation parameter, a, versus radius r (g+ is shown on the top, </_ 
in the bottom). Radius is in units of GM/c 2 . The horizon is shown in black. The curve a(r)t r — rme (dashed) 
is also plotted across the contour lines. Notice that, in the traditional disk-line scheme, no radiation is 
supposed to originate from radii r < r ms . If this is the case, then one must assume that all photons originate 
outside the dashed curve, and so the effect of frame-dragging is further reduced. 



butions from r > 10. For example, setting r = 7 
and O = 30°, one finds that the relative differ- 
ence between g + (a = 1) and g+(a — 0) amounts 
to only ~ 2%. A similar value comes out for g_, 
and this also roughly corresponds to the precision 
which is needed in order to be able to reveal the 
effect of black hole rotation. Notice, however, that 
this is in principle achievable by the Chandra grat- 
ings and the Astro-E2 calorimeter. It may also be 
achievable by XMM-Newton with high signal-to- 
noise data, and of course the effect would be more 
pronounced in case of a source seen edge-on. 

Other user-defined emissivities can be easily 
adopted. This can be achieved either by using 
the convolution component kyconv or by adding 
a new user-defined model to xspec. The latter 
method is more flexible and faster, and hence rec- 
ommended. In both approaches, the ray-tracing 
routine is linked and used for relativistic blurring. 
Naturally, by adding new emissivity laws one in- 
variably introduces additional freedom (for exam- 
ple, the height h of the primary source in the lamp- 
post model). It is worth repeating that caution is 
always necessary because the predictive power of 



the model rapidly decays with the number of pa- 
rameters. 

We conclude this section by summarizing the 
way of using the new model: Time-dependent 
spectra are computed by specifying the photon 
numbers emitted from the disk, N\ oc (E, r, ip, /xi oc , i), 
where r and ip are polar coordinates in the disk 
plane, /ii oc is cosine of the emission angle with 
respect to the disk normal, t is Boyer-Lindquist 
time in the Kerr metric, and N comes out in 
photons cm" 2 s" 1 sr _1 keV -1 . Naturally, mean 
spectra represent a special case in which time 
dependence of the source emissivity is not consid- 
ered. 

3. Application to MCG 6-30-15 

The Seyfert 1 galaxy MCG-6-30-15 is a unique 
source in which the evidence of a broad and skewed 
Fe Ka line has led to a wide acceptance of mod- 
els with an accreting black hole in the nucleus 
(Tanaka et al. 1995; Iwasawa et al. 1996; Nan- 
dra et al. 1997; Guainazzi et al. 1999; Fabian & 
Vaughan 2003). Being a nearby AGN (the galaxy 
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redshift is z — 0.0078), this source offers an un- 
precedented opportunity to explore directly the 
pattern of the accretion flow onto the central hole. 
The Fe K line shape and photon rcdshifts indi- 
cate that a large fraction of the emission originates 
from r < 10 (GM/c 2 ). The mean line profile de- 
rived from XMM-Newton observations is similar 
to the one observed previously using ASCA. The 
X-ray continuum shape in the hard spectral band 
was well determined from BeppoSAX data (e.g. 
Guainazzi et al. 1999). 

To illustrate the new ky model capabilities, 
we used our code to analyze the mean EPIC 
PN spectrum which we compiled from the long 
XMM-Newton 2001 campaign (e.g., as described 
in Fabian et al. 2002). The data were cleaned and 
reduced using standard data reduction routines, 
employing SAS version 5.4.I. 1 We summed the 
EPIC PN data from five observations made in the 
interval 11 July 2001 to 04 August 2001, obtaining 
a total good exposure time of <~ 290 ks (see Fabian 
et al. 2002 for details of the observations). The 
energy range was restricted to 3-10 keV with 339 
energy bins, unless otherwise stated, and models 
were fitted by minimizing the x 2 statistic. Statis- 
tical errors quoted correspond to 90% confidence 
for one interesting parameter (i.e. Ax 2 = 2.706), 
unless otherwise stated. 

No absorber was taken into account; the as- 
sumption here is that any curvature in the spec- 
trum above 3 keV is not due to absorption, and 
only due to the Fe K line. We remind the reader 
that the Fe K emission line dominates around 
the energy 6-7 keV, but it has been supposed to 
stretch down to ~ 3 keV or even further. We em- 
phasize that our aim here is to test the hypothesis 
that if all of the curvature is entirely due to the 
broad Fe K feature, is it possible to constrain a of 
the black hole? Obviously, if the answer to this is 
'no', it will also be negative if some of the spectral 
curvature between ~ 3-5 keV is due to processes 
other than the Fe K line emission. 

We considered Fe Ka and Fe K/3 iron lines 
(with their rest-frame energy fixed at E\ oc = 6.400 
and 7.056 keV, respectively), a narrow- line fea- 
ture (observed at E Q b s ~ 6.9 keV) modeled as 
a Gaussian with center energy, width, and in- 
tensity free, plus a power-law continuum. We 

1 See http://xmm.vilspa.esa.es/sas/. 



used a superposition of two kygline models to ac- 
count for the broad relativistic Ka and K/3 lines, 
zgauss for the narrow, high-energy line (which we 
find to be centered at 6.86 keV in the source 
frame, slightly redshifted relative to 6.966 keV, 
the rest-energy of Fe Lya), and kyhrefl for the 
Comptonized continuum convolved with the rel- 
ativistic kernel. We assumed that the high-energy 
narrow line is non-relativistic. Note that an al- 
ternative interpretation of the data (as pointed 
out by Fabian et al. 2002) is that there is not 
an emission line at ~ 6.9 keV, but He-like reso- 
nance absorption at ~ 6.7 keV. This would not 
affect our conclusions. The ratio of the inten- 
sities of the two relativistic lines was fixed at 
Ikp/Iku = 0.1133 (= 17/150), and the iron abun- 
dance for the Compton-reflection continuum was 
assumed to be three times the solar value (Fabian 
et al. 2002). We used the angles 6> min = 0° and 
#max = 90° for local illumination in kyhrefl (these 
values are equivalent to central illumination of an 
infinite disk). The normalization of the Compton- 
reflection continuum relative to the direct contin- 
uum is controlled by the effective reflection 'cov- 
ering factor', i? c , which is the ratio of the actual 
reflection normalization to that expected from the 
illumination of an infinite disk (see Appendix). 
Also r out (outer radius of the disk), a c (slope 
of power-law continuum radial emissivity in the 
kyhrefl model component), as well as R c , were in- 
cluded among the free parameters, but we found 
them to be only very poorly constrained. Weak 
constraints on a c and R c are actually expected in 
this model for two reasons. Firstly, in these data 
the continuum is indeed rather featureless. Sec- 
ondly, the model continuum was blurred with the 
relativistic kernel of kyhrefl, and so the dependence 
of the final spectrum on the exact form of emissiv- 
ity distribution over the disk must be quite weak. 

Model #1. Using a model with a plain power- 
law radial emissivity on the disk, we obtained the 
best-fitting values for the following set of parame- 
ters: a of the black hole, disk inclination angle Ol 
radii of the line-emitting region, r ln < r < r out , 
the corresponding radial emissivity power-law in- 
dex a of the line emission, as well as the radial 
extent of the continuum-emitting region, its pho- 
ton index T c , and the corresponding a c for the 
continuum. Notice that T c and a c refer to the 
continuum component before the relativistic ker- 
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Table 3: Spectral fitting results for MCG-6-30-15 using the ky model. 



# a 8 Continuum Broad Fe Ka line x 2 



£c etc Hn-Th r b -r h r out -r h q in q out EW (dof) 



1 0.35±o 3o 31.8 ±0.3 


2.01 ± 0.02 


1.0+.? 


5.1 ± 0.2 


11.4 ± 0.8 


3.9 ±0.6 258tl 6 3 


368.8 
(329) 


2a 0.99 ± 0.01 40.4 ± 0.6 


2.03 ± 0.02 


5.5l 6 2 


0.67 ± 0.04 3.35 ± 0.05 


4o±r 3 ° 


6.9±° l 9.7±° l 268 ± 13 


308.6 
(330) 


2b 0.72^0.30 28.5 ± 0.5 


2.01 ± 0.01 




0.65 ± 0.35 2.1 ± 0.2 f 
7.2 ± 0.2 


48+™° 


8 1 + 1 - 4 4 9+ 4 241 + 13 


313.5 
(330) 


2c 0.25 ± 0.03 27.6 ± 0.6 


1.97 ± 0.02 


"^-o.i 


1.23 ± 0.06 4 ±0.02 


109+.?°, 


9.2 ± 0.2 3.1 ± 0.1 267 ± 10 


313.9 
(330) 



Best-fitting values of the important parameters and their statistical errors for models #1-2, described in the text. The 
models include broad Fe Kct and Fe K/3 emission lines, a narrow Gaussian line at ~ 6.9 keV, and a Compton-rcflcction 
continuum from a rclativistic disk. These models illustrate different assumptions about intrinsic emissivity of the disk 
(the radial emissivity law need not to be a simple power law, but axial symmetry has been still imposed here). The 
inclination angle O is in degrees, relative to the rotation axis; radii are expressed as an offset from the horizon (in 
GM/c 2 ); the equivalent width, EW, is in electron volts. 

t Two values of the transition radius define the interval (»"b_,r b+ } where reflection is diminished. 



nel was applied to deduce the observed spectrum. 
As a result of the integration across the disk, the 
model weakly constrains these parameters. 

There are two ways to interpret this. A model 
which is over-parameterized is undesirable from 
the point of view of deriving unique model param- 
eters from modeling the data. However, another 
interpretation is that a model with a greater num- 
ber of parameters may more faithfully reflect the 
real physics and it is the actual physical situation 
which leads to degeneracy in the model parame- 
ters. The latter implies that some model parame- 
ters can never be constrained uniquely, regardless 
of the quality of the data. In practice one must ap- 
ply both interpretations and assess the approach 
case by case, taking into account the quality of the 
data, and which parameters can be constrained by 
the data and which cannot. If preliminary fitting 
shows that large changes in a parameter do not 
affect the fit, then that parameter can be fixed at 
some value obtained by invoking physically rea- 
sonable arguments pertaining to the situation. 

We performed various fits with the inner edge 
tied to the marginally stable orbit and also fits 
where r- m was allowed to vary independently. Free- 
fall motion with constant angular momentum was 
assumed below r ms , if the emitting region ex- 
tended that far. Table 3 gives best-fitting values 
of the key relativistic line and continuum model 
parameters for the case in which n n and a were 
independent. 



Next, we froze some of the parameters at their 
best-fit values and examined the x 2 space by vary- 
ing the remaining free parameters. That way we 
constructed joint confidence contours in the plane 
a versus O (see left panel of Fig. 3). These rep- 
resentative plots demonstrate that o appears to 
be tightly constrained, while a is allowed to vary 
over a large interval around the best fit, extending 
down to a = 0. 

Model #2. We explored the possibility that the 
broad-line emission does not conform to a unique 
power-law radial emissivity but that, instead, the 
line is produced in two concentric rings (a 'dual- 
ring' model). This case can be considered as a 
toy model for a more complex (non-power law) 
radial dependence of the line emission than the 
standard monotonic decline, which we represent 
here by allowing for different values of a in the 
inner and outer rings: a- m and a ou t- Effectively, 
large values of the power law index represent two 
separate rings. The two regions are matched at 
the transition radius, r = r b , and so this is essen- 
tially a broken power law. We explored both the 
case of continuous and discontinuous line emissiv- 
ity at r = 7V Notice that the double power-law 
emissivity arises naturally in the lamp-post model 
(Martocchia et al. 2000) in which the disk irradi- 
ation and the resulting Fe Ka reflection are sub- 
stantially anisotropic due to fast orbital motion 
in the inner ring. Although the lamp-post model 
is very simplified in several respects, namely the 
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Fig. 3. — Confidence contours around the best-fitting parameter values (indicated by a cross). Left: the 
case of a single line-emitting region (model #1 with zero emissivity for r < r ms ). Right: the case of a non- 
monotonic radial emissivity, model #2b. The joint two-parameter contour levels for a versus O correspond 
to 68%, 95% and 99% confidence. 



way in which the primary source is set up on the 
rotation axis, one can expect fairly similar irra- 
diation to arise from more sophisticated schemes 
of coronal flares distributed above the disk plane. 
Also, in order to provide a physical picture of the 
steep emissivity found in XMM-Newton data of 
MCG-6-30-15, Wilms ct al. (2001) invoked strong 
magnetic stresses acting in the innermost part of 
the system, assuming that they are able to dissi- 
pate a considerable amount of energy in the disk 
at very small radii. Intense self-irradiation of the 
inner disk may further contribute to the effect. 

This more complex model is consistent with the 
findings of Fabian et al. (2002). Indeed, the fit is 
improved relative to models with a simple emissiv- 
ity law because the enormous red wing and rela- 
tively sharp core of the line are better reproduced 
thanks to the contribution from a highly redshifted 
inner disk (see Table 3, model #2a). For the same 
reason that the more complex model reproduces 
the line core along with the red wing well, the 
model prefers higher values of a and O than what 
we found for the case #1. Notice that a — ► 1 
implies that all radiation is produced above r ms . 
Maximum rotation is favoured with both a and 
9 appearing to be tightly constrained near their 
best-fit values. Likewise for the continuum radial 
emissivity indices. There is a certain freedom in 
the parameter values that can be accomodated by 
this model. By scanning the remaining param- 
eters, we checked that the reduced x 2 <~ 1 can 



be achieved also for a going down to ~ 0.9 and 
8 <~ 37°. This conclusion is also consistent with 
the case for large a in Dabrowski et al. (1997); 
however we actually do not support the claim that 
the current data require a large value of a. As 
shown below, reasonable assumptions about in- 
trinsic emissivity can fit the data with small a 
equally well. 

In model #2a, small residuals remain near E <~ 
4.8 keV (at about the ~ 1% level), the origin of 
which cannot easily be clarified with the time- 
averaged data that we employ now. The excess 
is reminiscent of a Doppler horn typical of rela- 
tivistic line emission from a disk, so it may also be 
due to Fe Ka emission which is locally enhanced 
on some part of the disk. We were able to re- 
produce the peak by modifying the emissivity at 
the transition radius, where the broken power-law 
emissivity changes its slope (model #2b). We can 
even allow non-zero emissivity below r ms (the in- 
ner ring) with a gap of zero emissivity between 
the outer edge of the inner ring and the inner 
edge of the outer one. The inner ring, r; n <r<rb-, 
contributes to the red tail of the line while the 
outer ring, rb+<r<r out , forms the main body of 
the broad line. The resulting plot of joint confi- 
dence contours of o versus a is shown in Fig. 3 
(right panel) In order to construct the confidence 
contours we scanned a broad interval of parame- 
ters, < a < 1 and < 8 Q < 45°; here, detail is 
shown only around the minimum \ 2 region. 
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Two examples of the spectral profiles are shown 
in Fig. 4. It can be seen that the overall shapes 
are very similar and the changes concern mainly 
the red wing of the profile. For completeness we 
also fitted several modifications of the model #2 
and found that with this type of disk emissivity 
(i.e. one which does not decrease monotonically 
with radius) we can still achieve comparably good 
fits which have small values of a. The case #2c in 
Table 3 gives another example. This shows that 
one cannot draw firm conclusions about a (and 
some of the other model parameters) based on 
the redshifted part of the line using current data. 
The data used here represent the highest signal- 
to-noise relativistically broadened Fe K line profile 
yet available for any AGN. 

The differences in x 2 between model #1 and 
models #2a, 2b, and 2c are very large (~ 60) for 
only one additional free parameter, so all varia- 
tions of model 2 shown are better fits than the 
case #1. We note that formally, model #2a and 
#2c appear to constrain a much more tightly than 
model #1. The issue here is that the values of a 
can be completely different, with the statistical er- 
rors on a not overlapping (compare, for example, 
model #2a and #2c in Table 3). This demon- 
strates that, at least for the parameter a, it is not 
simply the statistical error only which determines 
whether a small a value is or is not allowed by 
the data. We confirmed this conclusion with more 
complicated models obtained in ky by relaxing the 
assumption of axial symmetry; additional degrees 
of freedom do not change the previous results. 

One should bear in mind a well-known techni- 
cal difficulty which is frequently encountered while 
scanning the parameter space of complex models 
and producing confidence contour plots similar to 
Fig. 3. That is, in a rich parameter space the pro- 
cedure may be caught in a local minimum which 
produces an acceptable statistical measure of the 
goodness of fit and appears to tightly constrain 
parameters near the best-fitting values. How- 
ever, manual searching revealed equally accept- 
able results in rather remote parts of the param- 
eter space. Indeed, as the results above show, 
we were able to find acceptable fits with the cen- 
tral black hole rotating either slowly or rapidly 
(in terms of the a parameter). This fact is not 
in contradiction with previous results (e.g. Fabian 
et al. 2002) because we assumed a different ra- 



dial profile of intrinsic emissivity, but it indicates 
intricacy of unambiguous determination of model 
parameters. We therefore need more observational 
constraints on realistic physical mechanisms to be 
able to fit complicated models to actual data with 
sufficient confidence (Ballantyne, Ross & Fabian 
2001; Nayakshin & Kazanas 2002; Rozariska et al. 
2002). 

We have seen that the models described above 
are able to constrain parameters with rather dif- 
ferent degrees of uncertainty. It turns out that the 
more complex type #2 models (i.e. those that have 
radial emissivity profiles which are non-monotonic 
or even have an appreciable contribution from 
t < r ms ) provide better fits to the data but a 
physical interpretation is not obvious. Ballantyne, 
Vaughan & Fabian (2003) also deduce a dual- 
reflector model from the same data and propose 
that the outer reflection is due to the disk being 
warped or flared with increasing radius. 

4. Discussion and Conclusions 

We have presented an extended computational 
scheme which can be conveniently employed to 
examine predicted X-ray spectral features from 
black-hole accretion disks. As mentioned above, 
a one-dimensional version of the code basically 
reproduces previous methods, with the addition 
of the option to vary more parameters (e.g. a). 
The current two-dimensional version also accomo- 
dates non-axisymmetric and time-evolving emis- 
sivity from a geometrically thin disk. The new tool 
facilitates accurate comparisons between model 
spectra and actual data. It provides also polari- 
metric information. One can even substitute other 
tables for the Kerr metric tables to explore dif- 
ferent specetime. However, we defer detailed dis- 
cussion of these latter capabilities to a future pa- 
per. Furthermore, work is in progress to include 
realistic emissivities relevant to different physical 
situations, in order to allow for the intrinsically 
three-dimensional geometry of the sources (e.g. an 
optically-thin corona above the disk), and to par- 
allelize the code. 

Recently, several research groups have em- 
barked on projects to compare theoretical pre- 
dictions of intrinsic emissivity of accretion disks 
with observational data obtained by X-ray satel- 
lites. This approach offers a fascinating opportu- 
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Fig. 4. — Spectrum and best-fitting model for the XMM-Newton data for MCG-6-30-15 in which the Fe K 
line originates in a dual-ring. The models 2a (left) and 2b (right) are shown in comparison. See Table 3 
for the model parameters. The continuum component (kyhrefl) is also plotted. The data points are not 
unfolded: the spectrum in these units was made by multiplying the ratio of measured counts to the counts 
predicted by the best-fitting model and then this ratio was multiplied by the best-fitting model and then by 
E 2 . 



nity to explore processes near black-hole sources, 
both AGNs and BHCs, and to estimate physical 
parameters of black holes themselves. In order 
to limit the number of unknown model param- 
eters and to alleviate the computational burden 
of the fit procedure, simplifications on the model 
emissivity (for example, stationarity and/or axial 
symmetry) have been often imposed. 2 Naturally, 
fitting data with complex models is computation- 
ally more demanding. 

To summarize, we have illustrated the code by 
following common practice and employing min- 
imisation of x 2 t° find the best-fitting parame- 
ter values for the highest signal-to-noise Fe K line 
profile available to date — from the XMM-Newton 
campaign of the Seyfert 1 galaxy MCG-6-30-15. 
The new model can vary more of the parameters 
than previous variants of the disk-line scheme al- 
lowed, which seems to be necessary in view of the 

2 Traditionally, Occam's razor is invoked to justify simpli- 
fying assumptions. However, this way of reasoning is not 
sufficiently substantiated here in view of the fact that ac- 
cretion flows are normally found to be very turbulent and 
fluctuating rather than smooth and steady. Put in other 
words, the assumption of an axisymmetric and stationary 
disk may provide a successful fit to the time averaged spec- 
trum, but it is not satisfactory because we know that more 
complicated models must be involved on physical grounds. 
("Things should be made as simple as possible — but not 
too simple. . . ") 



fact that simple prescriptions for intrinsic emisiv- 
ity (axisymmetric, radially monotonic, etc.) are 
not adequate to reflect realistic simulations of ac- 
cretion flows. It is without much surprise that a 
complex model exhibits an intricate \ 2 space in 
which ambiguities arise with multiple islands of 
acceptable parameters. The mean XMM-Newton 
spectrum of MCG-6-30-15 analyzed here justi- 
fies rejection of the often-used plain power-law 
emissivity that depends solely on radius from the 
black hole. However, we cannot arbitrate the 
issue of non-monotonic emissivity versus a non- 
axisymmetric one. Due to the inherent degen- 
eracy in time-averaged spectra, both alternatives 
can reproduce these data with reasonable parame- 
ters. In other words, a of the black hole cannot be 
unambiguously determined if complex (realistic) 
emissivities are adopted. 

It is clear that the relativistic line profiles af- 
fected by strong gravity are such that it is very dif- 
ficult to recover unambiguous information about 
the key parameters, such as black-hole spin, disk 
inclination angle, and inner disk radius. This is be- 
cause even with increased throughput, energy res- 
olution and signal-to-noise, the smooth, feature- 
less tails of the time-averaged line profiles do not 
retain sufficient information to break the degen- 
eracies between the dependences on the black- hole, 
accretion disk, and line emissivity parameters. 
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Higly specific spectral shapes could actually be 
provided by contributions to the line profile from 
higher-order images, whose form inherits unique 
information about space-time of the black hole. 
It may appear rather unlikely that higher-order 
images could be discovered in spectra of AGNs, 
because this would require to see an unobscured 
inner disk at large (almost edge-on) inclination. 

As the data improve, the important measurable 
quantities that will constrain parameters when 
compared with the models are: (i) the slope of 
the red wing, (ii) the energy of the cut-off of the 
blue side of the line profile, and (iii) the detailed 
shape of the blue peak (see Figs. 7-9 and further 
discussion in the Appendix). The extreme part of 
the red wing will not be as useful because it will 
always be difficult to tell where the line emission 
joins the continuum no matter how good the data 
are, and also because there is often no prominent 
sharp feature in the red wing since the red Dopplcr 
peak is too smeared out when the line emission 
extends too close to the black hole. In addition, 
time variability will be extremely important. No 
matter how, or if, the line emission responds to 
the continuum, we know that the disk inclination 
and black-hole spin very likely do not change dur- 
ing different snapshots. Therefore any differences 
in the snapshot line profiles must be due only to 
changes in the line emissivity, in terms of its ra- 
dial and azimuthal distribution and/or its radial 
extent (n n and r out ). Combining these approaches 
and collecting enough snapshots may provide suf- 
ficient information to uniquely determine one of 
the parameters, such as black-hole spin. 
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Newton data reduction and use of the Xassist 
software package. 3 We thank A. C. Fabian for 
helpful discussions and for providing us with a 
relativistic blurring routine which we used to 
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hole rotation, B. Czcrny for discussion about 
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the first version of the manuscript. A. Mar- 
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3 See http:/ /xassist. pha.jhu.edu/. 
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A. The Numerical Code: Tests and Examples 



0.00 0.05 0.09 0.14 0.18 0.23 0.00 0.20 0.40 0.60 0.80 1.00 




0.00 0.21 0.42 0.63 0.84 1.05 0.84 0.88 0.92 0.96 1.01 1.05 




Fig. 5. — Graphical representation of the contents of the data tables for the ky model used in the compu- 
tations. Several sets of contour plots capture the equatorial plane with an increasing resolution near the 
black-hole horizon. Here we show the effect of energy shift g=g(r,ip). Lengths are in GM/c 2 (radius up 
to r = 10 3 in the coarsest grid). Values are encoded by a colour scale, as indicated in each graph. Key 
parameters are: O — 30°, and a = 0.9987 (corresponding to rh=1.05, r ms =1.2). Clock-wise distortion of the 
contours is due to frame-dragging near a rapidly rotating Kerr black hole, and it is clearly visible in Boyer- 
Lindquist coordinates here. High accuracy of the tables has been achieved by eliminating this influence, with 
the aid of an appropriate transformation. In the upper left panel, corresponding to the region closest to the 
black hole, the horizon has been brought to the origin by the transformation of coordinates described in the 
text. 

Geometrical optics is adequate to formulate the ray tracing problem in rigorous manner (Schneider, Ehlcrs 
& Falco 1992). Physically, one has to consider a ray bundle connecting the source and the observer. We 
realize that the photon number is conserved in a bundle along the photon ray. Furthermore, this number is 
observer-independent quantity. Four-momentum of each photon is parallely transported along the ray. On 
a practical side, to achieve fast computations of observed spectra and for data fitting we produced accurate 
tables describing the Kerr metric, ray-tracing equations and the equations of geodesic deviation. These were 
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Fig. 6. — The same as in previous figure, but now showing the gravitational lcnsing effect in the disk. 
Contours correspond to lines of constant magnification. They were computed by solving geodesic deviation 
equations in the Kerr metric. 



solved in well-behaved coordinates. This choice allows one to map the whole region, from the outer edge 
of the disk down to the black hole horizon, where dragging effects are prominent. Also, caustics can play 
a role for large inclination (Rauch & Blandford 1994). When integration of light rays is accomplished, we 
transform results back to Boyer-Lindquist coordinates, so that the intrinsic emissivity of the disc can be 
conveniently defined as a function of r, tp and t in the disk plane. Here we briefly summarise the adopted 
approach. 

The source appears as a point-like object for a distant observer, so that the observer actually measures the 
flux entering the solid angle df2 b SJ which is associated with the detector area dS b s =D 2 dfi b s (this relation 
defines the distance D). We denote total photon flux received by a detector, 



^o 5 bs(£) 



dn(E) 
didS*, 



p-= fN loc (E/g)g 2 dn, 

obs J 



(Al) 
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Energy [keV] Energy [keV] 

Fig. 7. — Comparative examples of simple line profiles, showing a theoretical line (E\ oc = 6.4 keV) with 
relativistic effects originating from a black-hole accretion disk. Different sizes of the annular region (axially 
symmetric) have been considered, assuming that the intrinsic emissivity obeys a power law in the radial 
direction (a — 3). Resolution of the line-emitting region was n r xn v — 3000 x 1500 with a non-equidistant 
layout of the grid in Kerr ingoing coordinates, as described in the text. Left: dependence on the inner edge. 
Values of r; n are indicated in the plot (the outer edge has been fixed at the maximum radius covered by 
our tables, r out = 10 3 ). Right: dependence on r out (with the inner edge at horizon, r ln = r^). Other key 
parameters are: O — 45°, a = 1.0. Locally isotropic emission was assumed in the disk co-rotating frame. 



where 

AT (TP \ - dn l°c(^loc) , A9 ^ 

is local photon flux at the disk, dn(E) is the number of photons with energy in interval (E, E + dE) and 
g = E I E\ oc is the redshift factor. Note that Ni oc (E\ oc ) may vary over the disk, change in time, and it can also 
depend on the emission angle (for the sake of brevity of formulae, we do not always write this dependence 
explicitly). 

The emission arriving in solid angle dQ originates from area dS\ oc on the disk. Hence, in our computations 
we want to integrate the flux contributions over a fine mesh on the disk surface. To this aim, we adjust 
eq. (Al) to the form 



Here dS[ denotes an element of area on the disk perpendicular to light rays (corresponding to the solid angle 
dfl in flat space-time), and dS = rdrd^s stands for area in Boyer-Lindquist coordinates lying in the disk 
plane from which photons are emitted. The integrated flux per unit solid angle is 

dN? hs (E,t) ee NH s (E,t)dE 

/•r OVl t /-Vmax r(E+dE)/g 

= No dr dip dEi oc Ni oc {Ei oc ,r,ip,ni oc ,t Q )g 2 lr, (A4) 

J Tin t^min *> E / Q 

where to = t — At is time coordinate corrected for the light-time effect, 

dSi dSf 

l= dS^ c = dSl^ (A5) 

is the lensing factor in the limit D — > oo (keeping D 2 dfl constant), /ii oc is the cosine of the local emission 
angle, dS± is proper area at the disk perpendicular to light rays and corresponding to df2. In eq. (A4) we 
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Fig. 8. — More calculated line profiles, as in the previous figure. Left: Line profiles for different values of a. 
Notice the enhanced red tail of the line when the intrinsic emission is concentrated to the center of the disk. 
Right: Line emission originating from four different azimuthal segments of the disk. This plot can serve as a 
toy model of non-axisymmetric emissivity or to illustrate the expected effects of disk obscuration. Obviously, 
the receding segment of the disk contributes mainly to the low-energy tail of the line while the approaching 
segment constitutes the prominent high-energy peak. These two plots illustrate a mutual interplay between 
the effect of changing a and the impact of obscuration, which complicates interpretation of time-averaged 
spectra. The radial range is < r < 10 3 in both panels. 



used transformation of the proper area element, cL§i oc , 

= g~ L - (A6) 



dSioc _ i 



dS 

Finally, N^ s (E,t) is the output of the model computations. 

In eq. ( A4) ,No is a normalization constant of the model. For line emission, normalization is chosen in such 
a way that the total flux from the disk is unity. In the case of a continuum model, the flux is normalized to 
unity at certain value of observed energy (typically at E = 1 keV, as in other XSPEC models). Notice that 
default normalization has been changed in Figs. 7-9 for the sake of clarity. 

We employed the Bulirsch-Stoer method for numerical integration. Four sets of data tables were pre- 
computed and stored with the necessary information about photons received from different regions of the 
source. Specifically, separate tables hold information on (i) the energy shift (g-factor, or the ratio of energy 
of a photon received by an observer at infinity to the intrinsic energy when emitted at the disk), (ii) the 
lensing effect (the ratio I of the area subtended by photons at infinity, perpendicular to light rays through 
which photons arrive at the detector, to the area on the disk from where these photons originate), (iii) the 
mutual delay of the photons arriving at the observer's location (i.e. relative time lags At), and (iv) the local 
emission direction ^i oc with respect to the disk normal (evaluated in the disk co-rotating frame). 

The problem of directional distribution of the reflected radiation is quite complicated and the matter 
has not been completely settled yet (e.g. George & Fabian 1991; Ghisellini, Haardt & Matt 1994; Zycki 
& Czerny 1994; Magdziarz & Zdiarski 1995). A specific angular dependence is often assumed in models, 
such as the limb-darkening in laor xspec model, I(fJ,\ oc ) °c 1 + 2.06/j,i oc . However, it has been argued that 
limb-brightening may actually occur in the case of strong primary irradiation of the disk, and this is relevant 
for accretion disks near black holes where the effects of emission anisotropy are crucial (Czerny et al. 2004). 
We find that the spectrum of the inner disk turns out to be very sensitive to the adopted angular dependence 
of the emission, and so the possibility to modify this profile and examine the results using ky appears to be 
rather useful (we set locally isotropic emission as default). Impact of the darkening law on the reflection 
spectrum is connected with the assumed rotation law of the disk material. As mentioned above, we assumed 
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Fig. 9. — Dependence on observer inclination (values of 9 Q are given in the plot). Left: non-rotating black 
hole, a — 0. Right: maximally rotating black hole, a = 1. Other parameters as in Fig. 7. 



Keplerian circular motion for r>r ms . For r < r ms we assumed free-fall trajectories, maintaining constant 
orbital energy and angular momentum. 

All data tables are stored in fits format 4 and loaded in memory only once, when the routine is initiated. 
Examples of the contents of the data tables, illustrating the effects of frame-dragging, are shown in Figures 5- 
6 for given choices of a and O . The graphs represent a top view of the equatorial plane, each one on four 
different spatial scales. In practice we have chosen r out = 10 3 (Boyer-Lindquist radius) as a maximum outer 
limit of the disk, but the most dramatic dependences on, for example a, happen much closer to the horizon. 
In the upper left panels of Figs. 5 and 6 we define the radial coordinate in a different way to that in the 
other panels. Specifically, we use r' 2 = x' 2 + y' 2 = (r — r^) 2 . The use of r' brings the horizon to the origin 
so that the region just outside it is well resolved in the plots that show the region closest to the black hole. 
In these graphs, the observer is located at the top of the pictures with an inclination angle of 30° relative to 
the rotation axis. The marginally stable orbit is also shown (drawn as a circle) where relevant. 

Despite the fact that the particular choice of the graphical representation is a technicality to certain extent, 
we find these kinds of graphs useful for quick estimations of the expected range of energy shifts, fluctuations 
of radiation flux, time delay effects etc. For example, we notice that with the adopted (typical) value of 
o = 30°, the lensing cannot account for more than a few percent effect on observed count rates. Lensing, 
however, becomes important for edge-on orientation of the disk, as can be seen from the tables with larger 
values of O . We therefore arranged a systematic atlas of the tables for different values of angular momentum 
and observer inclination. Further, we produced their transformation from Boyer-Lindquist (r, tp, t) to Kerr 
ingoing coordinates (in which the effect of frame dragging is largely eliminated), and the same for other 
relevant quantities (time delay and the emission angle in the disk). It is worth noticing that different 
space-time metric and the disk rotation law can be accomodated simply by replacing the data tables. Also, 
transfer of Stokes parameters can be computed, and to this aim we computed data tables which are needed 
for polarimetry in strong gravity. However, this discussion remains beyond the scope of the present paper. 5 

Results of an elementary code test are shown in Figures 7-9. The intrinsic emissivity was assumed to 
be a narrow gaussian line (width cr=0.42FWHM — 5eV) with the amplitude decreasing oc r~ a in the local 
frame co-moving with the disk medium. No background continuum is included here, so these lines can be 
compared with similar pure disk-line profiles obtained in previous papers (e.g. Laor 1991; Kojima 1991) 
which also imposed the assumption of axially symmetric and steady emission from an irradiated thin disk. 
Again, the intrinsic width of the line is assumed to be much less than the effects of broadening due to bulk 

4 Flexible Image Transport System. See e.g. Hanisch et al. (2001) or http://fits.gsfc.nasa.gov/ for specifications. 
5 The corresponding graphs can be obtained from an internet site, http://astro.mff.cuni.cz/ky/. 
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Keplerian motion and the central gravitational field. Furthermore, Figures 10-11 compare model spectra of 
widely used xspec models. 

Typically, the slope of the low-energy wing is rather sensitive do the radial dependence of emissivity. 
Notice also that the laor model gives zero contribution at energy below O.lEioc = 0.64 keV and its grid 
has only 35 radial points distributed in the whole range 1.23<r<400. Therefore, in spite of very efficient 
interpolation and smoothing of the final spectrum, the laor model does not accurately reproduce the line 
originating from a narrow ring. Also, dependence on the limb darkening/brightening cannot be examined 
with this model, because the form of directionality of the intrinsic emission is hard-wired in the code, together 
with the position of the inner edge at r>r ms . 

The diskline model has been also used frequently in the context of spectral fitting, assuming a disk around 
a non-rotating black hole. This model is analytical, and so it has clear advantages in XSPEC. Notice, 
however, that the lensing effect is neglected. All this affects especially the spectrum near the black hole, 
where radiation is expected to be very anisotropic and flow lines non-circular. 
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Fig. 10. — Comparison of the output from XSPEC models for the disk-line problem: laor, diskline, and kygline 
(line 1 corresponds to the same limb-darkening law and a — 0.9982 as in laor; line 2 assumes locally isotropic 
emission and a — as in diskline). Left panel: o — 30°, r in — 6, r out = 7. Right panel: 9 = 70°, r in = 100, 
?~out = 200. Radial decay of intrinsic emissivity follows a = 3 power law. Default normalization of the model 
has been retained in this plot. 



B. Fast Model for Compton Reflection 

In this Appendix we describe the Compton-reflection model href I. This model has been a part of the 
xspec distribution for some time. We have used similar approximation for the Compton-reflection in the 
disk to produce the relativistic version, kyhrefl. 

The model has similar features to the well-known pexrav model (Magdziarz & Zdziarski 1995), which can 
often be prohibitively slow. There is a need for a fast reflection code for spectral fitting purposes, since the 
energy resolution of X-ray spectral data has been improved considerably in the last few years. Also, when 
convolving the reflection continuum with a kernel for relativistic energy shifting, computation speed is a 
critical factor. The price to pay is that href I makes use of a number of simplifying approximations. One of 
these is the elastic scattering approximation. Still, href I is accurate enough to be used for modelling data 
from observations of low-redshift sources with ASCA, Chandra, and XMM-Newton, where the source-frame 
spectrum does not extend beyond ~ 15 keV. (Notice, however, that this condition may be violated if the 
redshift of received photons is very large, e.g. near a black-hole horizon.) Another advantage of hrefl is that 
the model is completely analytic, and this is one of the major factors which makes it fast. 
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Fig. 11. — Comparison between the laor and kygline models. These plots are similar to previous figure, but in 
logarithmic scale and for three choices of the darkening law in kygline — (1) I oc 1 + 2.06/ii oc ; (2) / = const; 
(3) I oc log(l + l//ik>c)- Notice the impact on the calculated red-wing slope. Left panel: o — 30°; Right 
panel: 9 a = 70°. In both panels, r- m = r ms , r out = 400, a — 0.9982. 



Below we give an outline of the model, providing important details of the adopted approximations, and 
a direct comparison with the pexrav model. A general description of the approach can be found in Basko 
(1978) along with any standard radiative transfer textbooks (e.g. Chandrasekhar 1960). 




Fig. 12. — Left: Comparison between pexrav (solid) and hrefl*powerlaw (dashed) for two different inclination 
angles. Approximations are involved, making the model valid below E ~ 15 keV in the source rest-frame. 
Right: Comparison between the models kyhrefl (solid line) and hrefl*powerlaw (dashed line). Only the 
reflection component of the radiation flux due to the primary power-law continuum is plotted (in arbitrary 
units). Default model parameters were used with 9 = 30°. 

Let us consider a point source situated at a height h above the center of an optically-thick disk with inner 
and outer radii r; n and r out respectively. We define the cosines of the angles between the disk axis and lines 
connecting the inner and outer radii to the source as /x m i n and ^ max respectively, where 

^vCW and """'T^W <B1) 

The parameters actually required by the href I model are 9 min = arccos/i m j n and # max = arccos/i max . In the 
approximation of an infinitely large disk with no hole at the center, m i n — > and # max — > tt/2. 
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Now we adopt cylindrical coordinates with an origin coincident with the center of the disk face, and a set 
of rays emanating from the source and incident on the disk, making an angle arccos/xi = arccos [h/y/h 2 + r 2 ] 
with the disk normal at the point of interception. The incident specific intensity is 7;(r, ip). In the case of 
diffusive reflection from a plane-parallel atmosphere, the standard solution for the specific intensity of the 
reflected rays (as a function of energy, e) gives 

I{ ^ = M h ^B±mBtlim [keV cm" 2 s" 1 sr" 1 key 1 ]. (B2) 
4tt + /j, T 

Here, /z r is the cosine of the reflected rays with respect to the disk normal, A(e) is the single-scattering 
albedo (i.e. ratio of Thomson-scattering cross-section to the total scattering plus absorption cross-section), 
and H[n,\(e)] arc the Chandrasckhar iJ-functions for isotropic scattering. The angles are to be evaluated 
in the frame co-rotating with the disk medium, which is particularly important in the relativistic version of 
the code when the bulk orbital motion is very rapid close to the black hole. The assumption of isotropic 
scattering is one of the simplifications we employ: at photon energies much less than the electron rest-mass, 
the angular dependence of photon-electron scattering is not very strong. The albedo is defined as 

A(e) = 12 '1? Try ( B3 ) 

1.2 (Tt + er a b s (e) 

where ctt is the Thomson cross-section and cr a b s is the photoelectric absorption cross-section. The numerical 
factor 1.2 accounts for electrons from ionized helium (those from heavier elements are neglected due to their 
smaller numbers). 

Now we can reformulate equation (B2) in terms of actual fluxes. Let iVd(e) be the direct spectrum 
(photons cm~ 2 s _1 keV -1 ) of an X-ray source, observed at a large distance D. Further, di 7 ! denotes the 
energy flux within an energy interval de, which is incident on the face of the disk at radius r from the disk 
center: 

dF i = / i Mi= 4 f^t\ d % i [kcVcm-s-]. (B4) 
4-7r(r 2 + h z ) 

The reflected flux that is received by a detector having an effective area A from a disk surface element of 
area, r dr dtp, is then 

dF r = ii I I(fi I ,e)rdrdip 4tt ( — ) 4 [keV cm~ 2 s" 1 ]. (B5) 

\4:TtD z J A 

Here, 47r( J 4/47rZ) 2 ) is the solid angle subtended by the detector, r — /itan#i, and dr = hdOJ fif. 

Further, we denote the observed reflected spectrum N T (e) [photons cm~ 2 s _1 keV -1 ]. Since d_F r = 
ed7V r (e) ded^i r , we can substitute eq. (B4) into eq. (B2), then eq. (B2) into eq. (B5), and by integration 
over <p from to 2ir we find 

7V r (e) - iA(e) N d (e) f H^, A(e)] Hfa, A(e)]. (B6) 

For the if-functions we use an approximation (Basko 1978) which is good to ~ 8% accuracy in the range 
< A(e) < 1 and < n < 1: 

g|M,A(6)]~ jA=h • (B?) 

1 + v 3 C 1 — A ( £ )) m 

Next, we notice that for a given A(e), H[fii,X(ej\ does not vary much, so we can bring it in front of the 
integral, replacing the exact formula with the mean over the range fn = to \i\ = 1, 

H[X(e)} = [ ff[ M i,A(e)]d M i. (B8) 
Jo 
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For A(e) < 1 we obtain 



mm 



1 - y/l - A(e) 
\/3(l-A(e)) 



In l + y/3(l-X(e)) + 



1 



(B9) 



and for A(e) = 1 



H[X(e)} = 1 + § V3. 



(BIO) 



Finally, we integrate equation (B6) over m to find 



N T (e,fj, T ) 



\ N d (e,^) A(e)^ r In 



5[A(e)]J?[ Mr ,A(e)] 



(Bll) 



Mr ~l~ Mmax 



The function href I is encoded as a multiplicative component in xspec. When multiplied by the direct observed 
continuum, hrefl returns the reflected continuum, which is further multiplied by an effective reflection fraction, 
R c . The resulting contribution is added to the fraction R d of direct continuum flux. In other words, R c = 1 
corresponds to 2tt of the disk covering solid angle as seen by the X-ray source. Thus, 



where /x r has the meaning of observer inclination. The arguments of hrefl are the free parameters of the 
model (there are two more that are not listed explicitly here: they are the iron abundance relative to solar, 
and energy of the Fe-K edge). Notice that R c and R d are not independent, so no more than one of these 
two parameters is allowed to be free at any time of the fitting procedure. This way of coding the model is 
simply a matter of convenience for the user. 

Fig. 12 (left panel) compares the pure reflection spectrum produced by hrefl with the pexrav model using 
an input power law (photon index T = 2) for two inclination angles (/i r = 0.95 and /i r = 0.1, corresponding 
to o = 18.2° and 84.3° respectively; pexrav does not allow for values of fi r < 0.95). In hrefl the parameters 
6* m i n and 9 max were fixed at 0° and 90° respectively for the spectral fitting described in the present paper. 
Note that the relative iron abundance parameter for hrefl was set to 1.418 = 4.68 x 10~ 5 /3.3 x 10~ 5 because 
the two codes use different iron abundances. 

Apart from the adopted approximations, additional discrepancies arise between hrefl and pexrav due to 
the use of different photoelectric absorption cross-sections (hrefl uses older versions of cross-section values 
than pexrav and the difference is visible at low energies of the spectrum). Between 1-15 keV the errors are 
less than 20% (much better accuracy than this value is achieved for the facc-on inclinations). If the reflection 
continuum is combined with the direct continuum, for typical covering factors (R c ), the errors below the 
Fe-K edge are practically negligible. 

We conclude with the following concise summary of the approximations used in hrefl. The starting 
point of the model is the standard solution of the radiative transfer equation for reflection in a plane- 
parallel atmosphere. The simplifying assumptions are as follows, (i) Isotropic elastic scattering is assumed, 
(ii) Approximate analytical forms of the Chandrasckhar ff-functions are used, as described in Basko (1978), 
but the accuracy of the approximate formulae is better than <~ 8%. (iii) When integrating the incident 
continuum flux on the disk surface, the corresponding f/-function is taken outside the integral and replaced 
with its mean value, H . 

Fig. 12 (right panel) shows an example of the kyhrefl model spectrum, in which hrefl spectrum was 
smeared across the disk, kyhrefl can be interpreted as a Compton- reflection model for which the source of 
primary irradiation is near above the disk, in contrast to the lamp-post scheme with the source on axis. The 
approximations for Compton reflection used in hrefl (and therefore also in kyhrefl) are valid below ~ 15 keV 
in the source (and disk) rest-frame. 



7V d (e,^ r ) x href\(e m i n ,e maX7 9 ,R c ,R d ,e) ^ N d (e,^ r ) R d + R c 



Nd(e,Mr)_ 



(B12) 
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